<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/1999/REC-html401-19991224/loose.dtd">
<html lang="en">
<head>
	<meta http-equiv="Content-Type" content="text/html; charset=utf-8">
	<title>Simulations In An Axisymmetric Coordinate System Example (k-Wave)</title>
	<link rel="stylesheet" href="kwavehelpstyle.css" type="text/css">
	<meta name="description" content="Simulations In An Axisymmetric Coordinate System Example.">
</head>

<body><div class="content">

<h1>Simulations In An Axisymmetric Coordinate System Example</h1>

<p>This example provides a simple demonstration of using k-Wave for the simulation and detection of the pressure field generated by an initial pressure distribution within an axisymmetric coordinate system. It builds on the <a href="example_ivp_homogeneous_medium.html">Homogeneous Propagation Medium</a> and <a href="example_ivp_heterogeneous_medium.html">Heterogeneous Propagation Medium</a> examples.</p>

<p>
    <ul>
        <li><a href="matlab:edit([getkWavePath('examples') 'example_ivp_axisymmetric_simulation.m']);" target="_top">Open the file in the MATLAB Editor</a></li>
        <li><a href="matlab:run([getkWavePath('examples') 'example_ivp_axisymmetric_simulation']);" target="_top">Run the file in MATLAB</a></li>
    </ul>
</p>

<h2>Contents</h2>
<div>
	<ul>
        <li><a href="#heading2">Creating the k-space grid and defining the medium properties</a></li>
        <li><a href="#heading3">Defining the initial pressure distribution and sensor mask</a></li>
        <li><a href="#heading4">Running the simulation</a></li>
    </ul>
</div>

<a name="heading2"></a>
<h2>Creating the k-space grid and defining the medium properties</h2>

<p>Axisymmetric simulations are performed in an analogous fashion to those in two-dimensions. However, in an axisymmetric coordinate system, the x-dimension corresponds to the axial direction, while the y-dimension corresponds to the radial direction, as shown below. The coordinate system is rotationally symmetric about the x-axis, thus a point on the y-axis corresponds to a continuous circle in 3D space.</p>

<img vspace="5" hspace="5" src="images/kspaceFirstOrderAS_01.png" style="width:176px;height:250px;" alt="">

<p>Once the grid parameters have been defined, the medium discretisation is performed by <code><a href="kWaveGrid.html">kWaveGrid</a></code> in an identical manner to 2D simulations. However, in the radial dimension, the first grid point now corresponds to the grid origin, i.e., <code>y = 0</code>. In comparison, for <code><a href="kspaceFirstOrder2D.html">kspaceFirstOrder2D</a></code>, the Cartesian point <code>y = 0</code> is in the middle of the computational grid. For a heterogeneous acoustic propagation medium, the medium properties are given as two-dimensional matrices the same size as the computational grid. In this example, the properties are divided into two half-spaces.</p>

<pre class="codeinput">
<span class="comment">% create the computational grid</span>
Nx = 128;           <span class="comment">% number of grid points in the axial (x) direction</span>
Ny = 64;            <span class="comment">% number of grid points in the radial (y) direction</span>
dx = 0.1e-3;        <span class="comment">% grid point spacing in the axial (x) direction [m]</span>
dy = 0.1e-3;        <span class="comment">% grid point spacing in the radial (y) direction [m]</span>
kgrid = kWaveGrid(Nx, dx, Ny, dy);

<span class="comment">% define the properties of the propagation medium</span>
medium.sound_speed = 1500 * ones(Nx, Ny);       <span class="comment">% [m/s]</span>
medium.sound_speed(Nx/2:end, :) = 1800;         <span class="comment">% [m/s]</span>
medium.density = 1000 * ones(Nx, Ny);           <span class="comment">% [kg/m^3]</span>
medium.density(Nx/2:end, :) = 1200;             <span class="comment">% [kg/m^3]</span>
</pre>

<a name="heading3"></a>
<h2>Defining the initial pressure distribution and sensor mask</h2>

<p>As in two-dimensions, the initial pressure distribution is set using a matrix which contains the initial pressure values for each grid point within the computational domain. However, for the axisymmetric code, the initial pressure is rotationally symmetric about the first column in the matrix. For example, if the initial pressure is set to a single grid point in the first column (i.e., on the x-axis for y = 0), this corresponds to setting a point source within a 3D simulation. Conversely, if the initial pressure is set to a single grid point in any other column, this corresponds to setting a ring in 3D. In this example, the initial pressure is set in the shape of half a disc, which corresponds to a ball in 3D.</p>

<pre class="codeinput">
<span class="comment">% create initial pressure distribution in the shape of a disc - this is
% generated on a 2D grid that is doubled in size in the radial (y)
% direction, and then trimmed so that only half the disc is retained</span>
source.p0 = 10 * makeDisc(Nx, 2 * Ny, Nx/4 + 8, Ny + 1, 5);
source.p0 = source.p0(:, Ny + 1:end);
</pre>

<p>Again the sensor mask, which defines the locations where the pressure field is recorded at each time-step, can be given as a list of Cartesian coordinates, a binary mask, or the grid coordinates of two opposing corners of a rectangle. In this example, a Cartesian sensor mask in the shape of a circle is defined. A visualisation of the initial pressure and the sensor mask is shown below. Note, for the axisymmetric code, in the radial (y) direction, the perfectly matched layer (PML) is only applied to the outer edge of the domain. The default PML size is 20 grid points.</p>

<pre class="codeinput">
<span class="comment">% define a Cartesian sensor mask with points in the shape of a circle</span>
sensor.mask = makeCartCircle(40 * dx, 50);

<span class="comment">% remove points from sensor mask where y < 0</span>
sensor.mask(:, sensor.mask(2, :) < 0) = [];
</pre>

<img vspace="5" hspace="5" src="images/example_ivp_axisymmetric_simulation_01.png" style="width:560px;height:420px;" alt="">

<a name="heading4"></a>
<h2>Running the simulation</h2>

<p>The computation is started by passing the four input structures, <code>kgrid</code>, <code>medium</code>, <code>source</code>, and <code>sensor</code> to <code><a href="kspaceFirstOrderAS.html">kspaceFirstOrderAS</a></code>. By default, a visualisation of the propagating wave-field and a status bar are displayed. As the function runs, status updates and computational parameters are printed to the command line.</p>

<pre class="codeinput">
Running k-Wave simulation...
  start time: 26-Apr-2018 10:42:39
  reference sound speed: 1800m/s
  dt: 16.6667ns, t_end: 9.5333us, time steps: 573
  input grid size: 128 by 64 grid points (12.8 by 6.4mm)
  maximum supported frequency: 7.5MHz
  smoothing p0 distribution...
  calculating Delaunay triangulation...
  precomputation completed in 0.30384s
  starting time loop...
  estimated simulation time 3.14s...
  simulation completed in 3.6796s
  total computation time 4.043s
</pre>

<p>When the time loop has completed, the function returns the recorded time series at each of sensor points defined by <code>sensor_mask</code>. The ordering is again dependent on whether a Cartesian or binary sensor mask is used. A visualisation of the recorded sensor data is given below.</p>

<img vspace="5" hspace="5" src="images/example_ivp_axisymmetric_simulation_02.png" style="width:560px;height:420px;" alt="">

<img vspace="5" hspace="5" src="images/example_ivp_axisymmetric_simulation_03.png" style="width:560px;height:420px;" alt="">

</div></body></html>